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SUMMARY 


Approximate  maximum  likelihood  estimates  are  obtained  lor 
the  ordinates  of  a step  function  spectral  density  in  the  Gaussian 
case.  The  estimates  are  simply  integral  averages  of  the  periodogram 
over  the  frequency  bands  in  which  the  density  is  constant.  Whittle'. -> 
form  of  the  approximate  likelihood  is  used.  Results  are  given 
for  scalar  and  vector  processes. 


Some  key  words:  Stationary  Gaussian  process;  Step  function 

spectral  density;  Periodogram  average;  Approximate 
maximum  likelihood  estimates. 


Introduction 


1. 


Estimates  of  the  spectral  density  of  a stationary  process 
are  typically  constructed  by  smoothing  the  periodgram.  The  motiva- 
tion for  this  procedure  is  well-known.  For  a large  sample  size 
and  under  suitable  conditions  the  expectation  of  the  periodogram  is 
approximately  equal  to  the  spectral  density.  However,  the  variance 
does  not  tend  to  zero  as  the  sample  size  increases.  Averaging  of 
the  periodogram  reduces  the  variance  and  is  a method  which  can  be 
used  to  produce  a consistent  estimate.  At  the  same  time  the  averaging 
may  increase  the  bias. 

Spectral  estimates  constructed  by  periodogram  averaging  are 
essentially  estimates  obtained  by  the  method  of  moments.  They  are 
formed  by  replacing  population  covariances  by  their  sample  analogues 
and  then  smoothing  for  the  purpose  of  stabilizing  the  variance.  In 
this  paper  we  shall  assume  a Gaussian  model  and  show  that  periodogram 
averages  are  approximate  maximum  likelihood  estimates  of  certain 
spectral  parameters.  The  spectral  density  is  assumed  to  be  a positive 
step  function  with  a finite  number  of  distinct  ordinates.  These 
ordinates  are  the  spectral  parameters  to  be  estimated. 

To  derive  the  results  we  shall  use  the  approximate  Gaussian 
likelihood  proposed  by  Whittle  (1953  a, b,  195^)  for  the  analysis 
of  finite  parameter  linear  time  series  models.  The  use  of  a step 
function  spectral  density  in  Whittle's  approximate  likelihood  leads 
to  a simplification  which  permits  easy  exact  solution  of  the 
approximate  maximum  likelihood  derivative  equations.  The  resulting 


2. 


estimates  are  in  fact  estimates  of  the  average  amount  of  power  in 
each  of  a set  of  bands  determined  by  a finite  partition  of  the 
frequency  axis. 

The  approximate  likelihood  to  be  used  is  discussed  in 
*2  and  results  for  scalar  time  series  are  derived  in  §i>.  The 
vector  case  is  treated  in  1^. 


2.  An  Approximate  Likelihood  for  a Gaussian  Time  Series 

Let  {xt,  t=0,  +1,...}  be  a stationary  Gaussian  process 
with  mean  0 and  spectral  density  f ( X ) ( -n  < X < n ) . We  assume 
throughout  that  f(X)  is  positive.  The  problem  to  be  considered 
is  estimation  of  f(X)  or  of  parameters  associated  with  f ( X ) when 
an  observation  of  x = ( x^ , . . . , x,p)  ' is  available.  The  exact 
density  of  x is 

(27r)'^T|z|  exp(-®x,Z_1x),  (l) 


where  Z is  the  covariance  matrix  of  x. 

Let  [et,  t=0,  +1,...}  be  a sequence  of  independent  variables, 

p 

each  with  mean  0 and  variance  a . Then  the  process  {x^l  defined 

by 


xf  = Z 6 .et  , (t=0,  + 1,  . . . ),  (2) 

j=0  J 3 


where  = 1 and  z“_Q 
spectral  density  associated  with 


is  termed  a linear  process. 
(2)  is 


The 
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2 oo  ...  2 

f(X)  = — | Z 6 e1XJ|  ( -77  < X <7 t). 

2 77  j = 0 J 


It  is  common  to  assume  that  each  6^.  in  (2)  is  a function  of  a 
vector  of  parameters  about  which  inference  is  desired.  This  is 
the  case,  for  example,  when  {x^}  is  an  ARMA  (p,q)  process. 


q 

z 

k=0 
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ket-k 


( t=0,  + 1, . . . ), 


where  aQ  = 8 Q = 1. 


For  the  case  of  Gaussian  models  of  the  form  (2)  Whittle 
(1953a, 1954)  replaced  (1)  by 


( 2tto 


2)-ST 


exp 


fp  A^T 

4tt  ' -77 


Mil  4 

f(M  J 


= (2t r)-T  exp 


T f77 
477  ' -77 


lOg  f ( X ) 4- 


where 


(3) 


T 2 

l(X)  = -L-  | z x.elXt|  ( -77  < X < 77  ) 

2ttT  t=l  L 

is  the  periodogram. 

Two  considerations  are  involved  in  the  approximation  (3)- 

-1  * 

One  is  the  replacement  of  Z by  the  matrix  Z aelined  by 


« 

z 


_i_r  eu(»-t)f-i(x)a, 
(iTR  J-„ 


s,  t=l, . . . ,T 
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which  leads  directly  to  the  exponent  in  the  first  line  of  (3)- 

* 1 

It  can  be  shown  that  in  fact  Z - Z is  positive  semidefinite 
for  every  stationary  process  of  the  form  (2)  with  a positive  spectral 
density.  Moreover,  the  rank  of  Z - Z is  min{2  max(p,q),T}  for 
an  ARMA  (p,q)  process  and  is  T if  Cx t } of  the  form  (2)  is  not 
an  ARMA  (p,q)  process  with  p and  q both  finite.  (See  Shaman,  1976.) 

The  other  feature  of  the  approximation  (3)  is  the  replacement 
of  |Z|^  by  (o^)^T.  Note  that  the  second  line  of  (3)  follows 
from  the  first  by 


2 

o 


log{  27rf  ( X ) ) aX 


2t t -it 

P 

where  o is  the  variance  of  the  one- step-ahead  prediction  error 
and  is  the  same  as  the  variance  of  e^.  in  (2).  (See,  e.g.,  Grenander 
and  Rosenblatt,  1957,  p.  69. ) 

The  replacement  of  jz|^  by  (0^)*^  may  be  motivated  as 
follows.  It  is  well-known  (Grenander  and  Rosenblatt,  1957,  pp-  103-5) 
that  the  eigenvalues  of  the  Toeplitz  matrix  Z are  approximately 
given  by  equally-spaced  spectral  ordinates.  Then,  denoting  the 
eigenvalues  by  X^  (j=l,...,T),  we  can  write 


log  |£| 


iog  kj 


f 


and  the  latter  sum  is  approximately 


T r77 

2ir  ' -7r 


log{  2Trf(  X ) }dX  . 


5. 
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This  heuristic  argument  I'or  the  approximation  of  |£| 
is  not  dependent  upon  the  representation  of  in  (2),  but 

rTT 

holds  for  any  process  for  which  j ^ log  f(X)dA  > Moreover, 

_1 

the  same  argument  indicates  that  the  eigenvalues  of  £ are  ap- 

-1  * 

proximately  given  by  equally-spaced  values  of  f (X).  The  matrix  Z 
may  be  taken  as  an  approximation  to  Z for  processes  more 

general  than  (2).  In  fact,  Z~  = Z if  both  are  (two-sided) 
infinite-dimensional  matrices.  That  is, 

00  * 

Z o(s-r)o  (r-t)  = &(s-t)  ( s, t=0,+l, . . . ) , (4) 

r=  -oo 

where  &(•)  is  the  Kronecker  delta  and  o( • ) and  o (•)  are  the 

* 

elements  of  the  Toeplitz  matrices  £ and  Z , respectively.  The 
result  (4)  follows  because  the  left  side  of  (4)  is  the  convolution 
of  two  covariance  sequences,  and  the  Fourier  transform  of  this 
convolution  is  2tt  times  the  product  of  the  corresponding  spectral 
densities. 

For  the  above  reasons  we  shall  use  the  approximation  to  (1) 
given  by  the  second  line  of  (5)  when  {xt}  is  any  stationary  process 
with  a positive  spectral  density.  The  vector  case  analogue  of  this 
approximation  will  be  used  in  fi4. 

3-  Periodogram  Averages  as  Approximate  Maximum  Likelihood  Estimates. 

In  this  section  we  assume  the  spectral  density  is  a step 
function  and  derive  approximate  maximum  likelihood  estimates  for 


r 
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the  ordinates  of  the  density.  Let 

q 

J= 

where 


f(\)  = 2 f , c . ( X ) (0<X<tt), 

.i=l  J J 


(5) 


Cj(M  = l (^1  < * < V 

= 0 (otherwise). 


f ( -X ) = f(X),  and  0 = XQ  < \±  < . . . < Xq  = tt.  Set  f(rr)  = fq. 
Since 

r'V)  = x — o (X) 

j=i  t}  J 

and  l(X)  is  symmetric,  the  approximate  density  (3)  becomes 

r , ^ 

I n I «TT  *TT 

(2tt)”T  exp 


277  J 


q I pTr  1 rTr 

Z J log  f,  , c,(X)dX  + ! l(X)c,(X)dX\ 

=1)  J'0J  f “0  J 


(6) 


Let  L denote  the  corresponding  approximate  log  likelihood.  Then 


dL 
df . 


ptt  j c (X)  l(X)c  (X) 

1 1 -i- \ dX  (j=l,...,q). 


2n  0 


and  approximate  maximum  likelihood  estimates  are  given  by 

fTT 


fJ  - 


J0 

pTT 

j0  Cj(X)dX 


Wl 


X1 
r J 

jX  *(M<* 

j-1 
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( • • • > q)  • 
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Now  we  consider  asymptotic  properties  of  the  estimates  (7)- 
Care  is  required  in  the  analysis  because  t'(X)  is  not  continuous. 
First  we  note  that  <ff.  converges  to  f . as  T ->  ®.  To  see 

t]  u 

this  write 
X . 

<*!  ’ l(X)dX  - (X  - X )f  . 


X,  X . 

f Jr71'  p 0 

F,p(X-v)f (v)dVdX  - i f(X)dX 

’ xj-i  -w  Vi 


J Ft(v)  '■  tf(v  + X)  - f(X)}  dX  dv. 


where 


FT(v)  = 


. 21 
sin  a Tv 

rj 

2ttT  sin  ^v 


is  the  Fejer  kernel.  The  convergence  follows  as  in  Theorem  1.5.1 
of  Dym  and  McKean  (1972).  Although  f(X)  is  not  continuous,  the 
inner  integral  in  the  last  line  of  (8)  is  a continuous  function 
of  v.  From  this  fact  and  properties  of  the  Fejer  kernel,  we 
obtain  the  stronger  result  that  Ta(fj  - f j ) converges  to  0 for 
0 < a < 1. 

Next  consider 

T * (fj  - fj)2  = TVar  f + T(4fj  - fj)2. 


The  second  term  on  the  right-hand  side  tends  to  0 as  T-*®,  and 


/ 
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T(X1  -X1-i)  Var*i  = 


~Xj  X 


j j-i'  J *■  X . "X 


r J 

Cov{  I ( X ) , l( v ) } dA.dv 


j-l  j-l 


X,  X 
__  m f*  J f J 

K 

"j-l  J-l 


J T 

J -7T 


Ft(u  + X,M  -v)f(u)du 


(9) 


where 


dXdv, 


(x,v) 


sinfrXT  sin^vT 
2ttT  sin  JX  sin  sv 


(see  Anderson,  1971,  Theorem  8.2.8).  The  second  summand  in  the 

last  form  of  (9)  is 

Xj  X .-X  it  v 

t j J J J j Ft(x.x+  8)FT(y,y  + 8)r(x  - ^)f(y-*.)dxdy<l9(U. 

"j-l  "j-1'"  -n  -n 

A detailed  evaluation  of  the  limit  of  this  expression  is  lengthy  and 
tedious.  One  follows  the  same  steps  as  in,  e.g..  Theorem  9*3* 1 
or  Problem  16  of  Chapter  9 in  Anderson  (1971),  except  for  an  adjust- 
ment to  account  for  the  fact  that  f(X)  Is  not  assumed  to  be  con- 
tinuous here.  This  adjustment  consists  of  observing  that 

Xj 

J f(x-X)f(y-X)dX 
Xj-1 

is  uniformly  continuous  in  x and  y.  The  limit  of  the  first 
summand  in  the  last  form  of  (9)  is  0.  It  follows  that 


9- 


I 


lim  T(X  . - X * )2  Var  f = 2tt  f J f2(X)dX, 

T oo  J J X . , 

or 

ltm  T Var  f . 

T ^ oo  J 

Moreover, 

lim  T Cov  (f,,f.  ) = 0 (JA)> 

T °o  J 

We  make  several  comments.  First,  the  main  purpose  of  the 

calculations  in  this  section  has  been  to  demonstrate  that  periodogram 

averages  emerge  as  approximate  maximum  likelihood  estimates  of  the 

ordinates  of  a step  function  spectral  density  in  the  Gaussian  case. 

If  the  density  is  in  fact  a step  function,  it  is  unlikely,  though, 

that  the  statistician  desiring  to  estimate  the  ordinates  will  know 

the  correct  partition  (X  j=0,l, . . . , q) . For  frequencies  v and 

J 

u with  v>p  the  average 

1 rv 

— J l(X)dX 
v_u  p 

has  expected  value  which  converges  to  the  same  expression  with 
f ( X ) in  place  of  l(X),  and  T times  the  variance  converges  to 


However,  spectral  densities  encountered  in  practice  are  best 
viewed  as  being  continuous.  Thus,  the  model  (5)  really  should  be 
regarded  as  an  approximation  to  a continuous  density.  One  may 


2tt  f 


,2 


X .-X  . , 
J J-l 


( » • • • j q ) 
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interpret  this  approximation  as  having  arisen  from  smearing  of 
the  original  continuous  density.  The  model  (5)  will  be  difficult 
to  employ  in  practice  unless  a reasonable  partition  for  [ 0 , tt ] 
is  formulated  a priori. 

A final  comment  is  that  the  asymptotic  information  matrix 
corresponding  to  the  approximate  likelihood  (6)  is  diagonal  with 
elements 


Thus  the  estimates  f . are  asymptotically  efficient  in  this  sense. 

J 

4.  The  Vector  Case 

The  results  above  extend  to  the  vector  case.  Let  {x  ] be 
an  m-dimensional  Gaussian  process  with  mean  0 and  spectral 
density  £(x)  (-tt<x<tt).  The  matrix  £’(X)  is  Hermitian  and  is 
assumed  to  be  positive  definite  for  each  X. 

Denote  a time  series  by  ^ = ( x2  > • • • , x-p)  and 
x = vec(£)  be  the  mT  x 1 vector  obtained  by  arranging  the  columns 
of  J(,  from  left  to  right,  under  one  another  in  a single  column. 

If  £ denotes  the  mT  * mT  covariance  matrix  of  we  write 


11. 


/ s +•  \ 

The  mXm  submatrices  Zv  » ' are  approximated  by 


1 f 


( s , t- 1 , • • • »T) i 


(2ir) 


2 -TT 


' -1 


and  this  leads  to  approximation  of  z Z X by 


— J tr  X_1(X)i(X)dX, 

2tt  -it 


where 


i(X)  = 


T 

I e 


lX(s_t)  xsx^  (-Tr<A<m) 


2ttT  s,t=l 

is  the  periodogram.  The  covariance  determinant  |£|  is  ap- 
proximately the  determinant  of 


£(—),  £(—  ),■■■, X(2ir) 

T T 


2v  diag  ^ 

Thus  an  approximate  Gaussian  likelihood  for  £ is 


. -mT 

(2ir)  exp 


rv 


Ul  -TT 


log|f(X)|+  tr  f-1(X)l(X)l  dX 


(10) 


in  analogy  to  (3).  (See  Whittle,  1953b,  Theorem  6.)  Since  f(-X) 
= X(X)  and  I(-X)  = ^'(X),  we  may  rewrite  (10)  as 


mT 

exp 

r t rl 

J < 

2 log|£(X)|  +tr(f1(X)I(X)  +£"1(Ml'(*)) 

* 

i dX 

4tt  0 

m 

l J 

- 

:: 


. (ii) 
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Now  let 


f(\)  = Z Cj(X) 


J-l 


(0  < X < tt)  , 


(12) 


where  c.(X)  is  defined  below  (5)  and  the  Hermitian  matrices  f, 
J J 

are  positive  definite.  Denote 


~(tjk*  ) " ( fjki  + 1 fjk/ 


where  the  bar  designates  complex  conjugate  and  the  asterisk  conjugate 
transpose.  The  notation  {&)  designates  the  element  in  row  j, 
column  k of  the  matrix 

Substituting  (12)  into  (11),  we  obtain 


(27r)-mr  exp 


T 

Utt 


q 

z 

J=i 


2 log  |f.|  J c .(X)dX 

J UQ  J 


+ j tr^j1  l(X)  + c^.(X)dX 

If  L again  denotes  the  approximate  log  likelihood  and  is 

an  mXm  matrix  with  a 1 in  row  j,  column  k and  O' s elsewhere, 

then 


df 


dL 

jkk 


2w 


+ ~t  C trlfj1^£j1 1(k>  + fj'l£Kk£j'li,(X))  cj(x)dk 


- - — r cj<x)d>. 


x 

2? r 


+ — f XT1  KMX-1  + u:1  KMXi1}'  c (X)dX 

4tt  0 L J J J j kk  j 


( k— 1 9 • • • f m,  j— 


The  other  derivatives  are 


m A •R  P',T 

- (Xj  ) , j c.(X)dX 

TT  J k-*  0 J 


+ — J x:1  kmx:1  + U-1  kmx-1}' 

2?  0 J 3 J j 


c j ( X) dX, 


T / -1  T ^ 

1 (x,1)1  : 


; c (xjdx 
0 3 


"*  a 

T 1 1 4 >*  I 

- — j X1  KMxf  - fXiKMx;1)  C (X)dX 

2ir  0 J J J J kl  3 

m 

(Ml,  k,  1=1,  — ,m,  j=l, . . . , q) . 

“1  -1 

Since  l(X)fj  Is  Hermitian,  the  above  derivative  expressions 

may  be  simplified.  When  the  resulting  forms  are  set  equal  to  0, 
the  equations  become 

ij'.r  oj(»)d»  . J*  JijlejdXrt 

( J=i » • • • » q)  # 


f 


14. 


and  the  solution  L 


1.7/ 


o J 


r7r 

i 

0 


c .(^  )dX 


\ . 
r J 

J I(M^ 

Xj-i 


( j=l, . . . f q). 


The  asymptotic  properties  noted  for  the  scalar  case  generalize  to 
the  vector  case. 


References 

Anderson,  T.  W.  (1971).  The  Statistical  Analysis  of  Time  Series. 

New  York:  Wiley. 

Dym,  H.  and  McKean,  H.P.  (1972).  Fourier  Series  and  Integrals. 

New  York:  Academic  Press. 

Grenander,  U.  and  Rosenblatt,  M.  (1997).  Statistical  Analysis  of 
Stationary  Time  Series.  New  York:  Wiley. 

Shaman,  P.  (1976).  Approximations  for  stationary  covariance  matrices 
and  their  inverses  with  application  to  ARMA  models.  Ann. 
Statist.  4,  292-301. 


15- 

Whittle,  P.  (1953a)-  Estimation  and  information  in  stationary 


I 

time  eries.  Ark ■ Mat.  2_,  423-34. 

Whittle,  P.  (1953b)-  The  analysis  of  multiple  stationary  time 
serie:  . ■ J . R . Statist.  Soe.  B.  15.  125-39- 
Whittle,  P.  (1954).  Appendix  2 to  A Study  in  the  Analysis  of 
Stationary  Time  Series  by  H.  Wold.  Uppsala;  Almqvist  and 
Wiksell. 


' 


UNCLASSIFIED 


SECURITY  CLASSIFICATION  of  This  page  /When  CJete  Entered) 


REPORT  DOCUMENTATION  PAGE 


READ  INSTRUCTIONS 
BEFORE  COMPLETING  FORM 


l REPORT  NUMBER 

ONR  Report  No.  2 


2 GOVT  ACCESSION  NO 


3 RECIPIENT’S  CATALOG  NUMBER 


4 title  end  Submit) 


APPROXIMATE  I^XIMUM  LIKELIHOOD  ESTIMATION 
OF  A STEP  FUNCTION  SPECTRAL  DENSITY 


5 TYPE  OF  report  * PERIOO  COVERED 

7"  Technical  Rep^t. 


t performing  ORG.  REPORT  number 

Technical  Report  No.  11 


» CONTRACT  OR  ORRRT  NUMAERTi)' 

N00014-76-C-0950,  **Uj| 


7 Au  TmOR'IJ 

Paul /Shaman 


V 


9 PERFORMING  ORGANIZATION  name  and  address 

Department  of  Statistics 
Carnegie-Mellon  University 
Pittsburgh,  Pen n ylvania  15213 


T5  program  element,  project,  task 

, AREA  A WORK  UNIT  NUMBERS 

ACP'" 


/ 


NR-042-J567 


m 7 


II  CONTROLLING  office  name  and  address 

Office  of  Naval  Research 
Statistics  and  probability  Program 
Dock  (4jP).  Arlington.  Virginia  222 


fll 


Dii  REROAA 


Jun 


LL 


<1  NUMBER  OF  PAGES 

la 


* ^ /c(c 


11  SECURITY  cl-tst  (of  thla  rmpofi) 

Unclassified 


TJ  MONITORING  AGENCY  name  I AOORESSfll  dllterenl  troth  Controlling  Ottlco) 

T '/^'  US  OAJfcL 


IS*  OCCL  ASSlElCATlON  DOWNGRADING 
SCHEDULE 


>•  Distribution  STATEMENT  fo(  thla  Raport) 


Approved  for  public  release;  distribution  unlimited. 


17  Distribution  STATEMENT  fol  Iho  mbmttect  ontorod  In  Block  70,  II  dlttoront  Iroon  Report.) 


'•  Supplementary  notes 


It  KEY  BOROS  (Continue  on  rovoroo  oldo  II  nocoeoery  end  Identity  k y block  number) 


Stationary  Gaussian  process,  step  function  spectral  density, 
periodogram  average,  approximate  maximum  likelihood  estimates. 


To  AtilHACT  (Continut  of*  rovtao  oldo  If  ntcPltiR  aa*d  Identity  by  block  nmfe^r) 

pproximate  maximum  likelihood  estimates  are  obtained  for 
the  ordinates  of  a step  function  spectral  density  in  the  Gaussian 
case.  The  estimates  are  simply  integral  averages  of  the  periodo- 
gram over  the  frequency  bands  in  which  the  density  is  constant 

(continued  an  reverse  aide)  _ 
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Whittle's  form  of  the  approximate  likelihood  is  used, 
are  given  for  scalar  and  vector  processes.  JK 


Results 
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